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Introduction 


Modeling droplet condensation via CFD codes can be very tedious, time consuming, and 
inaccurate. CFD codes may be tedious and time consuming in terms of using Lagrangian 
particle tracking approaches or particle sizing bins. Also since many codes ignore 
conduction through the droplet and or the degradating effect of heat and mass transfer if 
noncondensible species are present, the solutions may be inaccurate. The modeling of a 
condensing spray chamber where the significant size of the water droplets and the time 
and distance these droplets take to fall, can make the effect of droplet conduction a 
physical factor that needs to be considered in the model. Furthermore the presence of 
even a relatively small amount of noncondensible has been shown to reduce the amount 
of condensation [Ref 1]. It is desirable then to create a modeling tool that addresses these 
issues. The path taken to create such a tool is illustrated. The application of this tool and 
subsequent results are based on the spray chamber in the Spacecraft Propulsion Research 
Facility (B2) located at NASA's Plum Brook Station that tested an RL-10 engine. The 
platform upon which the condensation physics is modeled is SINDA/FLUINT. The use 
of SINDA/FLUINT enables the ability to model various aspects of the entire testing 
facility, including the rocket exhaust duct flow and heat transfer to the exhaust duct wall. 
The ejector pumping system of the spray chamber is also easily implemented via 
SINDA/FLUINT. The goal is to create a transient one dimensional flow and heat transfer 
model beginning at the rocket, continuing through the condensing spray chamber, and 
finally ending with the ejector pumping system. However the model of the condensing 
spray chamber may be run independently of the rocket and ejector systems detail, with 
only appropriate mass flow boundary conditions placed at the entrance and exit of the 
condensing spray chamber model. The model of the condensing spray chamber takes 
into account droplet conduction as well as the degrading effect of mass and heat transfer 
due to the presence of noncondensibles. The one dimensional model of the condensing 
spray chamber makes no presupposition on the pressure profile within the chamber, 
allowing the implemented droplet physics of heat and mass transfer coupled to the 
SINDA/FLUINT solver to determine a transient pressure profile of the condensing spray 
chamber. Model results compare well to the RL-10 engine pressure test data. 

Facility and Exhaust System Description 

A brief description of the facility and its exhaust system is included in this section. 

Figure 1 is an aerial view of the Plum Brook Station B-2 test facility. 


4 



Figure 1: Aerial View of Spacecraft Propulsion Research Facility (B-2) 


Constructed in the 1960s, primarily to support the Centaur upper stage development, the 
Spacecraft Propulsion Research Facility (B-2) provides the facilities to simulate a space 
thermal soak and subsequent altitude firing of the propulsion system. Testing can consist 
of a variety of combinations including engine only, engine plus propellant delivery 
systems, or an integrated stage incorporating tanks and avionics. The facility is equipped 
with propellant delivery systems for LOX and LH2 plus helium and nitrogen supporting 
systems and is sized for hydrogen-oxygen engines up to 445 kN (100,000 lbf) thrust and 
approximately 200 kN (45,000 lbf) thrust for storable (non-condensable) propellant 
combinations. 

Space simulation is accomplished in a stainless steel cylindrical vacuum chamber 11.6 
meter (38 feet) diameter with a 18.9 meters (62 feet) vertical height. Vacuum pumping 
includes 3 stages of mechanical pumps and ten diffusion pumps ultimately bringing the 
vacuum chamber to a 10-4 Pa (10-5 Pa with liquid nitrogen in the cold wall) environment 
for well sealed systems. Thermal simulation is provided on the cold end by a liquid 
nitrogen cold wall and on the high end by portable lamps configured as required for the 
test2. 

Engine firing is accomplished by opening an 3.4 meter (11 ft) diameter valve [located at 
the end of the 12 meter (39 ft) diffuser] allowing the exhaust products to enter a spray 
chamber which cools and condenses the exhaust through circulation of 848 kL/min 
(224,000 gpm) of spray water from the water stored in the spray chamber. The 20.4 meter 
(67 ft) diameter by 36 meter (1 18 ft) deep concrete spray chamber is pumped by a steam 
ejector system to transport the remaining exhaust products to the atmosphere. 

To maintain vacuum conditions during engine firing, it is necessary to remove the 
products of combustion as fast as they are generated. At B-2, this is accomplished by a 
spray condensing system designed to condense steam (a product of combustion for 
hydrogen-oxygen engines) and an ejector system to pump out non-condensable gases 
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such as excess hydrogen. The ejectors are visible in the aerial view immediately to the 
right of the test building. 

Figure 2 is a cut-away view of the test building. A basic understanding can be made from 
this view as it shows a rocket upper stage (containing two engines) being lowered into the 
test chamber. Located immediately below the test chamber is the spray chamber where 
the rocket’s exhaust is exposed to water sprays (which condense a portion of the exhaust) 
as illustrated in Figure 3. Excess exhaust water vapor and the non-condensable gases are 
then removed by the primary ejector system which is powered by the steam supply 
system. The primary ejector system consists of 2 parallel trains with each train having 3- 
stages of ejectors. A quick observation reveals that much of the exhaust system’s 
function occurs below ground. 



Figure 2: Cutaway view of the B-2 test building 
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Figure 3: Condensing Spray System 


An important point to make is that as the test engines get larger and generate more 
products of combustion into the spray chamber, there will be some point where the spray 
chamber cannot condense enough of the rocket generated steam to maintain the low spray 
chamber pressures. The excess water vapor will start loading up the ejectors to the point 
where they cannot keep up. The spray chamber pressure will then begin a rapid pressure 
growth causing conditions at the engine to no longer satisfy test requirements. 


Model Description 

SINDA/FLUINT is used as the platform to model the rocket exhaust duct and spray 
chamber of the B2 facility. A schematic of the facility is shown in Figure 4. 
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The rocket exhaust is modeled using the Chemical Equilibrium with 
Applications code (CEA). CEA is included as a subroutine in the SINDA/FLUINT 
model and is therefore called from within the SINDA/FLUINT input file. The model is 
setup in a parametric manner so that all input, including CEA input, is set through the 
Register Data in SINDA/FLUINT. The input in the Register Data may also be input via 
an Excel Spreadsheet instead of through the ASCI SINDA/FLUINT input file. The 
expanded flow exiting the rocket, modeled in CEA, is input to the SINDA/FLUINT 
model of the rocket exhaust duct. The rocket exhaust duct flow is modeled using the 
FLUINT portion of SINDA/FLUINT (the fluids system analyzer) and the rocket duct 
wall is modeled in SIND A (the thermal analyzer). Thus heat transfer to the wall 
including any water spray used to cool the exterior of the duct wall may be modeled. The 
geometry of the rocket exhaust duct may be input in a simplified “lumped” approach or 
through more sophisticated methods, (i.e., full 3-D geometries). Currently the model has 
the ability to interface with MSC.Patran. The geometry can be built in MSC Patran, 
imported into SINDA/FLUINT, and then the thermal results can be mapped back on to 
the 3-D MSC Patran model. 
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Figure 5 though Figure 7 qualitatively illustrate the numerical model of the B2 facility as 
represented in SINDA/FLUINT. The model is composed of three primary submodels, 
“A”, “B”, “C”. Figure 5 illustrates Submodel “A”, which models the combined rocket 
duct flow and rocket duct wall model as described above. This submodel utilizes CEA to 
initialize the flow into the duct. Although Figure 5 illustrates fluid convection of the 
rocket duct flow and the rocket duct wall thermal conduction only in the radial direction, 
the wall geometry may include axial conduction as well. The MSC Patran model 
incorporates 3-D conduction through the rocket exhaust duct wall. 

At the end of the rocket exhaust duct CEA is called once more to calculate the shocked 
flow conditions. CEA is also used to model the quenching of the shocked flow. The 
quenched flow conditions are used as input to the spray chamber Submodel “B” of 
SINDA/FLUINT. The user may elect to run the SINDA/FLUINT model with only the 
rocket exhaust duct Submodel “A” and not with the spray chamber Submodel “B”. The 
shocking and quenching of the flow is contained in Submodel “A”, and updated in 
Submodel “B”. 

Figure 6 illustrates Submodel “B” which uses the shocked and quenched flow from 
Submodel “A” as input flow conditions to the spray chamber. The entrance location of 
the inflow is assumed to be at the bottom of the spray chamber which is also the location 
of the water line. The spray chamber is modeled in FLU1NT as a stratified vessel or a 
series of stacked “pancakes” that pass flow in one direction (upwards in Figure 56). The 
outflow of the spray chamber employs an ejector pump curve as a function tank pressure. 
The condensing spray flow dynamics and heat transfer is modeled via user FORTRAN 
and is based on particle flow mechanics. 

Figure 7 illustrates the thermal model of a single characteristic droplet that resides at a 
given time in a “pancake” as modeled in S1NDA. The heat rate on the particle is 
determined through user FORTRAN via a heat transfer coefficient which relates the 
thermal environment of its respective “pancake” determined by FLU1NT and the droplet 
wall temperature determined by S1NDA. The conduction through the droplet can be 
discretized to any number of “n” nodes. A radial conductor, “G n ”, may be determined via 
the following: 

Eq. 1 Gn = 4;ck d rnVl [W/K] 

- r . 
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Figure 7: SINDA Submodel “C” of Droplet 

It should be noted that Submodel “B” is inherently built as a true transient thereby 
making Submodel “C” also a true transient. This is because the “pancakes” illustrated in 
Figure 4Figure 6 are modeled in FLUINT as tanks which inherently have volume and 
therefore mass. The submodel’s physics require the use of these parameters. Submodel 
“A” however runs with the steady state solver but can be run as a “false” transient with a 
series of steady state runs. 

Model Details: Submodel “A” 

The SINDA/FLUINT Submodel “A” models the rocket duct flow and heat transfer 
through the rocket duct wall. The submodel also utilizes CEA to initialize the rocket duct 
flow, shock the exit duct flow, and quench the exit duct flow with cooling water. CEA is 
run within SINDA/FLUINT as a subroutine. Other applications of CEA will be outlined 
below as well as in the section, CEA Modeling Applications. 

SINDA/FLUINT Modeling of Supersonic Flow 

The rocket exhaust duct flow or duct entrance flow is supersonic. In order for 
SINDA/FLUINT to model supersonic flow several modeling techniques are employed. 
There are five significant issues that need to be addressed. These issues regard the use of 
FLUINT. First, a FLUINT set mass flow rate connector (MFRSET), is placed at the duct 
exit. Second, all choking calculations must be turned off in FLUINT. Third, set IPDC=0 
for the FLUINT connectors, i.e., duct friction calculations are supplied by the user. This 
is necessary for several reasons. Fluid properties used in the calculation of the friction 
factor must be evaluated at a reference temperature because of the extreme temperature 
variations which are produced in the compressible boundary layer. Currently, FLUINT 
does not evaluate fluid properties at a reference temperature in calculating friction 
factors. Equation 2 is used to evaluate fluid properties in calculating the friction factor 
[Ref 2]: 


11 



Tref- 0.5(Twall + Tfluid) + 0.22(T rec - Tstat) 


[K] 


Eq. 2 

The recovery temperature is defined by [Ref 3]: 

Eq. 3 Tree = Pr'^ (Tstag- Tstat) + Tstat [K] 

In order for FLUINT to obtain a stable supersonic solution, the friction factor must be 
input as a positive “FC” value with an “FPOW” exponent equal to one. Normally the 
“FC” value is negative. A negative value causes FLUINT’s solver to pursue a subsonic 
solution. The friction factor for flow through a duct with fluid properties evaluated at T re f 
is given by [Ref 3]: 


Eq. 4 F = 0.184Re~ a2 ^ 


The FLUINT “FC” factor can then be updated through user logic via the equation [Ref 
4]: 


Eq. 5 FC = — — [1/m/kg] 

2Ac"p 

It is important to note that in using this procedure to obtain a supersonic solution, that the 
DUCT MACRO option in FLUINT should not be used. This is because this macro 
automatically imposes an AC factor which is not appropriate for this case. Equivalent 
results can also be obtained by imposing a negative FK factor, the negative of Equation 4, 
or imposing a positive HC factor (pressure gain), Equation 5 multiplied by the mass flow 
rate squared. 

The fourth issue regarding supersonic flow involves heat transfer. The heat transfer to 
the duct wall is also supplied via user logic. A turbulent heat transfer coefficient is 
calculated with fluid properties evaluated at T re f using the Colburn Analogy [Ref 3]: 

Eq. 6 h D = 0.23 Re 08 Pr^ — [W/nf/K] 

E^d 

Since the convecting fluid temperature to the duct wall should be the recovery 
temperature, Eq. 2, and not the static temperature the following modeling techniques are 
employed. Thermal boundary nodes that analogously represent the fluid lumps in the 
duct are created. They are updated in user logic to correspond to the respective recovery 
temperature of the fluid lump. Thermal conductors are created from these thermal 
boundary nodes to the duct wall nodes. These conductors are updated in user logic using 
Eq. 6. The heat rate leaving a fluid lump and going into the respective duct wall node can 
be calculated using Eq 6: 
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Qd llD Asd (Tree — Twall ) 


[W] 


Eq. 7 

It should be noted that FLUINT does have a utility that incorporates a reference 
temperature in heat transfer calculations. The above procedure just illustrates an 
alternative approach. 

Finally, the fifth issue regarding FLUINT and supersonic flow is not trivial. Currently 
FLUINT imposes a velocity limit on the kinetic energy term in the total enthalpy energy 
equation for stability reasons. This limit does not allow for the maximum velocity 
attained in the supersonic flow modeled in the analysis of the B2 facility. The mach 
numbers calculated in the analysis of the B2 facility ranged from 0.0 to 8.0. For the exact 
velocity limit used in FLUINT, Cullimore & Ring Technologies should be consulted for 
this value may depend on the version of SINDA/FLUINT being used. The FLUINT 
maximum velocity constraint in this analysis was 3000 m/s. This constraint did not allow 
for the conservation of total enthalpy for adiabatic flow. In order to “conserve” total 
enthalpy a heat rate was imposed on the fluid lumps representing the duct flow. This heat 
rate represented the “pseudo” kinetic energy term that was missing because of the 
velocity limit. The heat rates for the individual fluid lumps were calculated in user logic 
so that the total enthalpies for the respective lumps were conserved. This process was 
also repeated for flow that included heat transfer to the wall. The entire process for 
determining the “pseudo” heat rates involved several SINDA/FLUINT FASTIC calls 
(i.e., steady state calls). First a steady state solution is obtained for the duct flow with no 
wall heat transfer and no imposed “pseudo” kinetic energy heat rates. This yields a 
solution that does not conserve total enthalpy. A second steady state call is then made 
that employs user logic to calculate the necessary heat rates to conserve total enthalpy for 
the individual fluid lumps. This second run restarts from the results of the first run. If 
wall heat transfer is modeled, a third steady state call is made that restarts from the results 
of the second steady state run. This third steady state call recalculates new “pseudo” heat 
rates that take into account the heat rate to the wall. 

In order to execute the aforementioned procedure, the stagnation conditions of the duct 
flow fluid lumps must be calculated. For no wall heat transfer this calculation need only 
be executed once since stagnation enthalpy is conserved. However for the case of wall 
heat transfer each fluid lump will have a new stagnation state. The problem in 
calculating stagnation conditions lies in that most property tables do not go to such high 
temperature conditions. One way to circumvent this is to use CEA to calculate stagnation 
conditions. 


CEA Modeling Applications 

CEA, Chemical Equilibrium with Applications, is a NASA developed code that 
calculates mixture chemical equilibrium compositions and properties. The source code is 
written in ANSI standard FORTRAN, and is appended as a subroutine to the 
SINDA/FLUINT model of the B2 facility. CEA is used for several functions in the 
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model. First, CEA determines the rocket exhaust duct flow properties that serve as inlet 
conditions to the rocket duct SINDA/FLUINT model. The rocket setup in CEA is run as 
an enthalpy/pressure case. For a given area ratio, CEA outputs parameters that can be 
used to calculate the rocket duct mass flow rate. The mass flow rate, pressure, and 
temperature output from CEA, for a given area ratio, are used as input to the 
SINDA/FLUINT rocket duct model. This is accomplished via user logic. 

CEA is also used to determine duct flow stagnation properties. For a given fluid lump 
where pressure, temperature, enthalpy, and entropy are known via SINDA/FLUINT, 

CEA may be used to determine the stagnation conditions of the fluid lump. A subroutine 
was created that sets up a CEA run and is called from user logic. The user may elect to 
send this subroutine the enthalpy, entropy, and pressure of a given fluid lump, determined 
via SINDA/FLUINT, or the user may elect to send the temperature and pressure of a 
given fluid lump, determined via SINDA/FLUINT. If temperature and pressure values 
are selected to be sent to the subroutine, then the subroutine first calls a 
temperature/pressure case in CEA so that the enthalpy and entropy may be determined. 
Also required as input to this subroutine, is the flow velocity associated with the fluid 
lump. Given the velocity and pressure of the fluid lump, a “guessed” value of the total 
pressure is calculated. The velocity and enthalpy, either supplied by SINDA/FLUINT or 
calculated via CEA, is used to calculate the total enthalpy of the fluid lump. The value of 
total pressure is used to run an entropy/pressure case in CEA. The output enthalpy 
calculated via CEA is then compared to the total enthalpy calculated from the user inputs 
of enthalpy and velocity. If these two values are not converged, a new guess on the total 
pressure is calculated and the procedure is repeated. Upon convergence, the subroutine 
returns total conditions of the fluid lump, (i.e., temperature, pressure and enthalpy). 

Another subroutine was created that sets up a CEA shock run. This subroutine is called 
from user logic in SINDA/FLUINT, using the fluid properties of the exit fluid lump in 
the rocket duct. For a given fluid lump where pressure and temperature are known via 
SINDA/FLUINT, as well as the velocity associated with the fluid lump, CEA calculates 
the post-shock flow properties. 

The properties of the shocked flow are used in the quench calculations. This shocked 
rocket duct flow exhaust is quenched with cooling water that runs along the outside of the 
rocket duct. Another subroutine utilizing CEA was created to calculate the quench water 
flow rate needed to bring the shocked rocket duct flow exhaust to the inlet quench 
cooling water temperature at the spray chamber pressure. Ultimately this is a transient 
process since the spray chamber pressure varies as a function of time. The total enthalpy 
of the shocked rocket exhaust duct flow, calculated using the aforementioned subroutine, 
as well as the enthalpy of the quenched cooling water, supplied via SINDA/FLUINT, is 
sent to this subroutine. It should be noted that the subcooling of the quenched cooling 
water is taken into account in the SINDA/FLUINT water property database. This 
subroutine sets up an iterative enthalpy/pressure calling procedure to CEA that guesses an 
initial mass flow rate of the quench cooling water. The enthalpy that is sent to CEA in 
this enthalpy/pressure procedure is defined by: 
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Eq. 8 


H 


quench 


ITlducl H (k|ct + mwaterH water 

• • 

lHwater lUduct 


[J/kg] 


The mass flow rate of the quench cooling water is iterated upon until the temperature of 
the vapor mixture, calculated by CEA, converges to the inlet quench cooling temperature. 
If the mass flow rate of the quench cooling water is less than the total inlet quench 
cooling water mass flow rate, the excess water is removed from the system (i.e., the water 
pool, see Figure 4). This quenched vapor serves as the inlet conditions to the spray 
chamber model, Submodel “B”. 

A final note regarding CEA and SINDA/FLUINT needs to be addressed. Many of the 
SINDA/FLUINT property databases are generated with REFPROP, Reference Fluid 
Thermodynamic and Transport Properties. This is a NIST standard reference database 
[Ref 5]. Care needs to be taken in transferring enthalpy and entropy values between CEA 
[Ref 6], and SINDA/FFUINT because the reference states used by each code may not 
coincide. Please refer to Appendix A for more information. 

Duct Wall Cooling 


The outside of the rocket duct is cooled by the quench cooling water. The quench 
cooling water is maintained at the spray chamber pressure and the inlet condenser spray 
temperature. A heat transfer coefficient is calculated that includes film boiling and 
convection. The film boiling heat transfer coefficient, which takes into account 
subcooling of the water, is given by [Ref 7]: 


Eq. 9 


where 


film 


kygPy (P£ -Py)hf g 

v Pvkk\i UC [ (To^uct T sat j ^ 


Eq. 10 


h[ g =(1.0 + 0.68 Ja)h fg 


and Ja is the Jakob number. 


[W/nr/K] 


[J/kg/K] 


The forced convection heat transfer coefficient is calculated via a SINDA/FFUINT 
subroutine, EXTCYF, which calculates the forced convection heat transfer coefficient for 
a semi-infinite circular cylinder subject to a perpendicular external flow [Ref 3]. In user 
logic the total heat transfer coefficient is calculated as [Ref 3]: 

Eq- 1 1 h quench = h film + h conv [ W/m 2 /K] 


Any water vapor formed by the quench cooling water may be optionally added to the 
spray chamber model, Submodel “B”. Placing this additional vapor above the spray bar 
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in the spray chamber will create a greater load on the ejector pumps and cause the 
pressure in the spray chamber to rise. 

Transient Control 


For a transient model of the rocket exhaust duct flow the series of FASTIC steady state 
calls, as previously described (SINDA/FLUINT Modeling of Supersonic Flow) , are 
called as desired by the user. A time increment for updating any rocket transient data is 
input through the SINDA/FLUINT Register Data Block. This transient data may include, 
but is not limited to, oxidizer to fuel weight ratio, and chamber pressure. The data may 
be input through the Array Data Block. The shocking and subsequent quenching of the 
rocket exhaust duct flow is respectively updated as well. However, the shocking and 
quenching of the duct flow may also be updated via a separate parameter during the 
transient analysis. The user may update this process as a function of change in spray 
chamber pressure via user logic. This incremental change of spray chamber pressure is 
input in the Register Data Block. 


Model Details: Submodel “B” 

Model Setup 


The B2 spray chamber is modeled in FLUINT as a stratified vessel using a series of 
stacked “pancakes”, previously shown in Figure 5. The “pancakes” are modeled as 
FLUINT tanks with ifaces between them. The tanks have volume and therefore mass, 
and the ifaces are “massless” membranes that will expand or contract the ta nk s to which 
they are attached, as mass is either added or subtracted from said tanks. These ifaces 
maintain constant pressure between the tanks. Also between the tanks are set volume 
flow rate connectors (VFRSET). Logic was created to calculate the appropriate volume 
flows rate to pass between the “pancakes” so that the volumes of the “pancakes” would 
remain relatively equal, or more realistically, within a certain percentage of each other. 
These volume flows rate always flow “up” and never reverse, that is, flow always 
proceeds in the direction from the rocket exhaust duct plenum to the ejector pump 
plenum. 

The shocked and quenched flow of the rocket exhaust duct enters the bottom of the spray 
chamber. In this analysis, this exhaust duct flow is comprised of steam and hydrogen. 
This inflow condition is modeled via a plenum with a set mass flow rate connector 
(MFRSET). These conditions are allowed to vary according to transient input data, as 
described in the previous section, Model Details: Submodel “A”, Transient Control . 
The outflow of the system is controlled by an ejector pumping system. This is modeled 
with a plenum and MFRSET connector. The ejector pumping capacity is input as mass 
flow rate versus pressure curves in the Register Data Block. These curves are used in the 
logic blocks to determine the appropriate amount of steam and hydrogen to remove from 
the system based on the spray chamber pressure. The ejector MFRSET connector utilizes 
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the species specific suction option to remove the appropriate amounts of steam and 
hydrogen based on the calculations using the ejector curves in the logic blocks. 

The model contains lumps which are “dummy plenums” that don’t have any impact on 
the flow solution but are used for bookkeeping purposes. These lumps store the saturated 
pressure and the core temperature of the characteristic liquid droplet that occupies a 
respective tank at a given time. Every ta nk has associated with it a “dummy plenum” that 
represents a droplet with a unique inlet Sauter Mean Diameter (SMD). For example if 
there are five stratified “pancakes”, or tanks, and only one inlet SMD family of droplets, 
then there are five “dummy plenums”. If there are two families of droplets with different 
inlet SMDs then there are ten “dummy plenums”, or two “dummy plenums” for every 
tank. Figure 8 illustrates the B2 spray chamber modeled as a stratified vessel using five 
FFUINT tanks. There is only one inlet SMD family of droplets so there are only five 
“dummy plenums” which represent the properties of the characteristic droplet within the 
respective five tanks. Although referred to as plenums, these FFUINT lumps are input as 
ta nk s but are held in boundary states via FFUINT subroutine calls to HFDFMP. This 
procedure essentially makes the tanks behave as plenums. This is done once again for 
bookkeeping purposes. Before the calls to HFDFMP are made, the volumes of these 
tanks are updated to store the current liquid volume of the sum total of the droplets 
occupying it respective stratified tank. 

The spray bar is modeled via a plenum with a set mass flow rate connector (MFRSET). 
Note that the location of the MFRSET connector will vary depending on how many tanks 
with which the user opts to stratify the spray chamber vessel. That is, the spray bar 
MFRSET will connect to the droplet “dummy plenum” whose respective stratified tank is 
at the approximate physical location of the spray bar. The spray bar plenum and its 
MFRSET connector, serve only to book keep the inlet water flow that enters the spray 
chamber. 
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Figure 8: FLUINT Model Setup of Spray Chamber 
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Figure 9: FLUINT Model Setup of Spray Chamber 

Each tank and droplet “dummy plenum” has associated with it two plenums, refer to 
Figure 9, which illustrates this for only one tank and its respective droplet “dummy 
plenum”. These plenums are updated in user logic to be at the saturation pressure and 
temperature of the characteristic droplet that occupies that tank at a given time. One 
plenum represents the saturated liquid state, and the other plenum represents the saturated 
vapor state. Note that this saturation state may not necessarily be equivalent to the 
saturation state of the tank. The saturated liquid plenum has an MFRSET connector 
which removes the condensate mass flow rate from its respective tank, and an MFRSET 
connector which places the condensate mass flow rate into its corresponding droplet 
“dummy plenum”. The latter MFRSET is just for bookkeeping purposes and does not 
influence the solution. The saturated vapor plenum has an MFRSET connector which 
transfers the evaporative mass flow rate from the droplets to its respective tank, and an 
MFRSET connector which removes the evaporative mass flow rate from the droplets 
from its corresponding droplet “dummy plenum”. The latter MFRSET is just for 
bookkeeping purposes and does not influence the solution. The saturated vapor plenum’s 
MFRSET connector may optionally include vapor formed from the quench cooling water 
used on the outside wall of the rocket exhaust duct. The total vapor flow rate formed 
from the quench cooling water may be equally divided and distributed amongst all the 
tanks, or it may be distributed equally amongst the tanks that physically reside above the 
spray bars. 

The shocked and quenched rocket exhaust duct flow rate that enters the spray chamber 
and passes through the stratified tanks never technically mixes with the spray bar inlet 
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water mass flow rate on the FLUINT level. That is, the spray bar water flow never 
occupies a FLUINT ta nk that contains the exhaust duct flow. The exhaust duct flow and 
the spray bar water flow follow counter flow paths that virtually run parallel to one 
another (see Figure 8). The water flow from the spray bar always flows “down” and is 
not modeled using MFRSET connectors because the water flow is not modeled as a 
continuum, but as droplets. The modeling of the droplet physics, momentum and heat 
transfer, is primarily accomplished through user logic. This is because it is desirable and 
more appropriate for the system being analyzed that the droplet mechanics be Lagrangian 
based. The FLUINT “pancake” is homogeneous and may be in thermal nonequilibrium 
with the droplets that reside in it at a given time. It is desired that a model be created that 
determines the condensing efficiency and not require this parameter as an input. Neither 
should the pressure distribution through the spray chamber ever be assumed but 
calculated in the model. Although the modeling of the droplet is Lagrangian based, 
tracking of individual droplets is not done. Rather, within each ta nk a time weighted 
average of droplet parameters is determined. These time weighted parameters represent 
the characteristic droplet that resides in a given tank at a given time. Also, for a given 
ta nk there is a characteristic droplet for each inlet SMD family of droplets. 

Droplet Model: Momentum 

The inflow water flow rate from the spray bar is modeled as an inflow of water droplets 
at a specified velocity. The droplets enter the spray chamber and interact with the 
FLUINT tank that physically resides at the approximate height of the spray bar. The 
exhaust duct flow through the VFRSET connectors does not include the drag force 
exerted upon it from the water droplets. However, the momentum of the water droplets 
does include the drag force exerted upon the water droplets from the exhaust duct flow. 
The relative velocity of a droplet that virtually resides in a FLUINT tank or “pancake” is 
defined by: 

Eq. 12 y rel =V d +V 0O [m/s] 

where V® is the velocity of the exhaust duct flow entering the VFRSET connector of the 
respective FLUINT tank (see Figure 10). 
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(Flow Upward) 


Positive Direction f 
Negative Direction | 


Figure 10: Characteristic Droplet in FLUINT Stratified Tank or “Pancake” 


Voo 




Figure 11: Forces Acting on a Droplet 

The forces acting on the droplet are illustrated in 

Figure 1 1 . The body force is defined by: 

Eq. 13 FG = m d g [N] 

or 

Eq- 14 FG= ^7ir 3 (p QO -p d )g [N] 
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The drag force is determined via a drag coefficient defined by: 


Eq. 15 


CD = 


-2.0FD 
^ r d 2 PooV rel |v rel | 


The drag coefficient is determined via Schiller Naumans relation [Ref 8] which is 
modified to have the correct limit for the inertial regime: 


Eq. 16 


rn 24.0 

CD = max 

Re 


((>■ 


0 + 0. 15 Re 


0.687 


),0.44) 


Figure 1 1 illustrate that the model always has the rocket exhaust duct flow through the 
VFRSET connectors go “up” or the positive direction, and the droplets flow “down” or 
the negative direction. From the FLUINT perspective, this is to avoid excessive flow 
reversals through the VFRSET connectors. User logic forces flow “upward” through a 
VFRSET connector if the volume of a given “pancake” is larger than the corresponding 
“pancake” above it, otherwise the flow rate through the VFRSET connector is zero. The 
droplets also do not experience flow reversal because of the complexity of book keeping 
them. If a characteristic droplet residing in a “pancake” at a given time experiences a net 
upward force the absolute velocity of the droplet will decrease until it reaches a value of 
zero. Once it reaches zero and continues to have a net upward force, the velocity 
continues to remain zero until there is a net downward force. During this time the 
characteristic droplet remains in its “pancake”. This in effect models “floating” droplets. 
Droplets from a “pancake” above a “pancake” with floating droplets are still allowed to 
enter the “pancake” where droplets are floating, given that the droplets in the above 
“pancake” are experiencing a net downward force. This will ultimately create an 
accumulation of droplets in the “pancake” still experiencing a net upward force with 
floating droplets thus causing that “pancake” to flood. The net force on the droplet is the 
sum of the body force and the drag force. The change in velocity of a droplet for a given 
time step is calculated via: 


Eq. 17 FG + FD = m d — [N] 

dt 

When a characteristic droplet travels the height of a “pancake” (see Figure 10) in which it 
currently resides at a given time, it will virtually pass to the “pancake” below. This 
means that all the droplets in the “pancake” move en mass to the “pancake” below. 
Because of this effect, there may be moments in time that a given “pancake” may be void 
of droplets. Another consequence of the droplets only being able to move “downward” is 
that droplets cannot virtually reside in any “pancake” above the spray bar. 

The number of droplets entering the spray chamber through the spray bar is calculated 
knowing the Sauter Mean Diameters (SMD) of the inlet families of droplets. At a given 
time the total number of droplets that virtually reside in the “pancake” that is 
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approximately at the physical height of the spray bar can be determined. Each new set of 
droplets that enter this “pancake” are at the inlet spray bar conditions and have a uniform 
temperature. The droplets that already reside in this “pancake” may have a different 
velocity and a different temperature distribution. Since each individual droplet is not 
tracked, a characteristic droplet having a time averaged value of velocity and a time 
averaged temperature distribution must be determined. The characteristic droplet has a 
velocity determined via: 


Eq. 18 


Vin d Ndrop in + (Vold d +dV d )Ndrop 
Ndrop in +Ndrop 


[m/s] 


where dVd is detennined from Equation 17. The temperature distribution will be 
discussed in Model Detail “Submodel B”, Droplet Model: Heat Transfer . 

Droplet Model: Heat Transfer 


The heat transfer between the droplets in a “pancake” and the exhaust duct flow involves: 

• sensible heat (conduction and convection) through the diffusion layer to the edge 
of the condensate layer over the droplet 

• the latent heat released from the vapor of the exhaust duct flow condensing over 
the droplet 

• sensible heat (conduction and convection) through the condensate layer over the 
droplet 

• thermal conduction through the droplet (see Figure 4) 

Figure 4 illustrates these heat transfer mechanisms. 



Figure 12: Heat Transfer Mechanisms between Exhaust Duct Flow and Droplet 
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The sensible heat (conduction transfer and convection transport) is modeled via Ranz 
Marshall [Ref 1]: 

Eq. 19 Nu =2. 0 + 0. 6 Re 0 ' 5 Pr 033 

where Nu is the Nusselt Number. The heat transfer coefficient is determined by: 

Eq. 20 Nu = h conv D d 

k 

The fluid properties are taken at a reference temperature using a simple 1/3 rule: 

Eq.21 T ref =T satd + tkfkal [K] 

The latent heat released from the water vapor in the exhaust duct flow condensing on the 
surface of a water droplet is determined by how much of this vapor is transported through 
the diffusion barrier of the noncondensible gas, hydrogen, building up next to the 
condensate layer. Reviewing Fick’s Law, [Ref 9], for a binary system of hydrogen and 
water vapor in spherical coordinates the transfer of water vapor is: 

Eq- 22 j w =~P D wh — ^ [kg/s/m 2 ] 

or 

and the transfer of hydrogen is: 

Eq- 23 j h = - p D wh — ^ [lcg/s/nr] 

or 

The mass flux rates of water vapor and hydrogen are respectively: 


Eq. 24 mw=Pw V w [kg/s/m 2 ] 

Eq. 25 m h =a V h [kg/s/m 2 ] 

The total average mass flux is: 

Eq. 26 pV=Pw V w+a V h [kg/s/m 2 ] 

and the mass averaged mixture velocity is: 

Eq. 27 V = x w V w +x h V h [m/s] 

The mass flux rates of water vapor and hydrogen consists of a diffusive part which is due 
to the relative velocity of the respective constituent to the mass averaged velocity of the 
mixture, and a convective part which is due to the “drift” of the respective constituent 
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with the mass averaged velocity of the mixture. Thus the absolute mass flux rates of 
water vapor and hydrogen respectively become: 


Eq. 28 

" 5x " " 

m w - pD wh + x w ( m + m w ) 

or h 

[kg/s/m 2 ] 

Eq. 29 

_ 5x. . " ", 

m h - pD wh + x h (m h + m w ) 

[kg/s/m 2 ] 

The above equations 

are subject to the condition: 


Eq. 30 

X w+ X h =1 



Under steady state conditions the mass flow rate of water vapor in the radial direction 
does not change: 


Eq. 31 


5(m w 47tr 2 ) = 8(m w ) 
Sr Sr 


[kg/s/m] 


The noncondensable gas, hydrogen, is taken to be insoluble in liquid water (i.e., at the 
surface of the water droplet). Under steady state conditions this translates to the mass 
flux rate of hydrogen being zero: 


Eq. 32 m =0 [kg/s/m 2 ] 

h 

If Equation 32 is imposed upon Equation 29 the result is: 

Sx 

Eq- 33 P D wh -^= x h m w [kg/s/m 2 ] 

In other words, during condensation, a driving diffusive force of hydrogen away from the 
droplet surface must exactly counterbalance the convective flow of the hydrogen toward 
the surface [Ref 10]. During condensation the noncondensable hydrogen gas 
accumulates at the surface of the water droplet. The partial pressure of the hydrogen at 
the surface of the droplet actually becomes larger than the partial pressure of the 
hydrogen in the bulk or exhaust duct flow. This build up of noncondensable gas has a 
significant effect on the heat transfer because the mass transfer of the water vapor is 
decreased because of the diffusion barrier. Also since the total pressure remains constant, 
at the droplet surface where the hydrogen partial pressure has increased the partial 
pressure of the water vapor falls below the saturation pressure of the water vapor in the 
bulk or exhaust duct flow. In turn this reduces the saturation temperature, T satd , at which 

condensation occurs. 

Applying Equation 32 to Equation 28 results in: 
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Eq. 34 


[kg/s/m 2 ] 


-pD 


wh 


ITLw 


8x w 

8r 


1 — x„ 


Applying Equation 3 1 to Applying Equation 34 results in: 

-47rr 2 pD wh 8x ” 

Eq. 35 m w = — 


8r 


1 — x„ 


[kg/s] 


or 


Eq. 36 


m 


A -4 Tip D.dx 


1 — x„ 


[kg/s/m] 


Integrating from the droplet surface to infinity, i.e. the exhaust duct flow, we obtain 


Eq. 37 


m* 

r 


= 4 K P D wh ln(l-x w )|~ 


[kg/s/m] 


or 


Eq. 38 


m, 


= 4 ^ r sP D wh ln 


"l-X A 

A WOO 

v‘- x wiy 


[kg/s] 


The mass flux rate then becomes: 


Eq. 39 


m v 


pD^in 


1 — x,. 


1-xi 


[kg/s/m 2 ] 


A negative value for the mass flow rate or the mass flux rate corresponds to “suction” or 
condensation on the droplet surface. The density in Equation 39 is evaluated at the 
reference temperature defined in Equation 21. If a Lewis number of one is assumed, i.e. 
the ratio of thermal diffusivity to mass diffusivity is equal to one, Equation 39 may be 
replaced by: 


Eq. 40 



k 

Cpr s 


In 


1 — x. 


1 — X wi J 


[kg/s/m 2 ] 


The heat transfer rate from the “bulk” or rocket exhaust duct flow to the droplet interface 
is: 
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Eq. 41 



[W] 


In a quasi-steady process this heat rate must equal the heat rate through the condensate 
layer. In the condensate layer the rate of heat transfer may be defined as: 


In order to solve forT satd , Equation 41 and Equation 42 must be iteratively solved. 
Basically one would guess an initial T satd , usually a value between T sat00 andT walld . Once 
T sa td ' s l <nown the mass fraction x wi may be determined. 

The condensation heat transfer coefficient is that for condensation on a rigid sphere [Ref 
1] with the following assumptions: 

• the temperature distribution in the condensate film is linear 

• axial conduction is neglected 

• shear forces at the interface are neglected 

Chung and Ayyaswamy, [Ref 1], make use of a correlation that includes inertial effects, 
(usually only film conduction is considered): 


Within the SINDA/FLUINT model an iterative technique similar to the one described in 
the previous section to calculate T satd is accomplished using a SINDA/FLUINT 

subroutine, HTUDIF, [Ref 4]. This subroutine uses the general purpose Chilton- 
Coulbum analogy: 


Eq. 42 



[W] 


Eq. 43 



[W/nr/K] 


Droplet Model: Heat Transfer Applied in SINDA/FLUINT 


Eq. 44 


h, 


conv 
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The interfacial temperature, T satd , is iteratively solved for using the following energy 
balance: 



For a given “pancake” or stratified ta nk the subroutine requires the following input: 

• Temperature of the mixture in the “pancake” 

• Absolute pressure of mixture in the “pancake” 

• Wall temperature (i.e., T walld ) 

• Mixture heat transfer coefficient (i.e., h conv HTUDIF ) 

• Condensation heat transfer coefficient (i.e., h cond ) 

The subroutine returns h gff , the effective heat transfer coefficient between the vapor 

saturation temperature in the “pancake” and the droplet wall temperature. Equation 45 is 
then used to calculate T satd . HTUDIF is called within the SINDA/FLUINT logic block, 

FLOGIC 0. It is important to note that Submodel “B” runs as a full transient. The 
implementation of the aforementioned logic implies a quasi-steady solution of the droplet 
interfacial temperature (droplet saturation temperature) at a given time step. 

The amount of condensate formed in a “pancake” at a given time is determined from 
Equation 43: 


This mass flow rate is removed from its respective “pancake” via an MFRSET connector 
with a species specific option to remove only water vapor. This MFRSET is illustrated in 
Figure 9 as the condensation MFRSET between Tank 2001 and Plenum 11001. However 
the MFRSET is removing water vapor at the absolute temperature of Tank 2001 and not 
at the actual interfacial temperature (droplet saturation temperature) T satd . The 

consequence of this is that the correct mass rate is being removed from the tank but not 
the correct energy rate and therefore a correction must be made to the energy input to the 
tank. The energy that is being removed from said tank via the MFRSET connector is: 


Eq. 46 



fg 


[kg/s] 


Eq. 47 


CU, = mcond H woo Ndrop 


[W] 


It should be: 


Eq. 48 


Q. c = mcond H satv Ndrop 


[W] 
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A heat rate must thus be applied to the tank. This applied heat rate should be: 


• • 

Eq. 49 QL = (m cond H woo - mcond H satv ) Ndrop [W] 

Basically Equation 49 “replaces” the energy lost in the ta nk through the MFRSET 
connector and removes the correct energy from the tank. 

Finally the convective heat rate leaving the tank must be included. The correct heat rate 
on a “pancake” or stratified ta nk is: 

• • 

Eq. 50 QL= [mcond H woo - mcond H satv - As s h COIlv (T o0 -T satd )]Ndrop [W] 

It should be noted that the amount of condensate formed is equally distributed among the 
number of droplets. The net increase in droplet diameter is negligible but book kept 
nonetheless. 

An issue that arises in the model occurs when the solution reaches steady state 
conditions. Implementation of the aforementioned logic is done at the beginning of a 
time step, and does not change as a solution is iterated upon within the internal solver of 
SINDA/FLUINT. The condensation rate that occurs within a “pancake” is constant and 
based on all the droplets that reside in that “pancake” as the solver finds a solution for 
that time step. The heat and mass transfer rates do not reflect the change in the partial 
pressure of the water vapor of the “pancake” that occurs as the solver attempts to find a 
solution for the given heat and mass transfer rates assigned at the beginning of the time 
step. Thus the partial pressure of the water vapor in the “pancake” may drop below the 
droplet saturation temperature because of the inability to correct for the heat and mass 
transfer rates as the partial pressure of the water vapor changes. The overestimation of 
the heat and mass transfer rates must be accomplished on the next time step, at which 
condensation will not occur but evaporation will. The mass flow rate of evaporation is 
similar to that of Equation 38: 


Eq. 51 


mevap — As s CF 


Cpr s 


-In 


1-P„ 


1-P ■ 

V W1 7 ) 


[kg/s] 


However, the mass fractions are replaced by partial pressures and a correction factor was 
added: 


Eq. 52 


CF=2.0 + 0.6Re Pr 


0.5 „ 0.33 


This correction factor accounts for the enhancement of evaporation due to convection 
[Ref 11]. 
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Instead of flow leaving the “pancake”, water vapor flows into the “pancake” via an 
MFRSET connector. This MFRSET is illustrated in Figure 9 as the evaporation 
MFRSET between Ta nk 2001 and Plenum 12001. Plenum 12001 is updated every time 
step to be the saturated vapor state of the characteristic droplet residing in Tank 2001. 

The heat rate on the “pancake” now becomes: 

Eq. 53 QL=[-As s h conv (X, - T atd )] Ndrop [ W] 

Also since condensation occurs at T satd which is lower than the saturation temperature of 
the water vapor in the FLUINT “pancake”, a very small amount of liquid may fonn when 
over-condensation occurs. If the quality becomes less than one, andT sat is less thanT satd 
the heat rate on the “pancake”, is modified to be: 

Eq- 54 QL = [— As s h conv (T^ - T satd ) + As s h film ^ (T satd - T sat ) (1 - XL)] Ndrop [W] 
where 


Eq. 55 




= 0.62 


I4gp v (ft-p v )hf g 

Pv D d( T satd- T sat) 


, 0.25 


[W/nr/K] 


The second term in Equation 54 accounts for the vaporization of any liquid mass fraction 
that may occur in the “pancake”. 


As the solution approaches a steady state the spray chamber pressure, as will be seen in 
the Results section, will tend to oscillate about a small pressure dead band. These 
oscillations can be minimized by reducing the FLUINT control variable, DTSIZF, 
thereby reducing the amount a lump property can change in a time step. Another way to 
minimize the oscillations is to use one of SINDA/FLUINT’s more robust matrix solvers, 
i.e. setting MATMET = 12. Finally one can increase the number of “pancakes” or 
stratified tanks. The aforementioned procedures may create longer run times and smaller 
time steps. 

Model Details: Submodel “C” 


Model Setup 


It is important to accommodate the physics of thermal conduction through the droplet in 
the SINDA/FLUINT model of the B2 facility. It is hypothesized that given the droplet 
sizes (on the order of 1500 microns and greater), droplet velocities (on the order of 37 
m/s), and size of the spray chamber, that the water droplets may not be fully utilized. 
Therefore SINDA “solid conduction” droplets that correspond to each of the time 
averaged characteristic droplets are modeled. A SINDA droplet that models the thermal 
conduction through the droplet corresponds to, in time and space, the characteristic 
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droplet of a FLUINT “pancake” discussed in the section, Model Details: Submodel 
“B”, Momentum . 

The number of thermal droplets created in the model equals the number of “pancakes” in 
the FLUINT model multiplied by the number droplet families (initial droplet sizes with 
varying SMDs). The degree of nodalization of the thermal droplet is determined by the 
user. For the current analysis the thermal droplets have ten nodes. The conductors 
between the nodes are defined by Equation 1 . All conductors and capacitances are 
updated as a function of time and temperature in user logic. 

These thermal droplets track the characteristic droplets created in the user logic blocks of 
the FLUINT portion of the model. That is, the time averaged properties of the droplet are 
overlaid onto the thermal droplet on the beginning of a time step in the thennal model. 
For the thermal model, the property of interest is the temperature distribution through the 
droplet. At the end of the time step the thermal model has solved for a new temperature 
distribution through a given droplet. This is the new temperature distribution for the 
characteristic droplet that is sent to the FLUINT user logic blocks. In the FLUINT logic 
blocks the droplet momentum equations are executed to determine whether droplets have 
left a given “pancake” and entered the “pancake” below. It is after this sequence of logic 
that a time average of droplet properties occurs. Equation 18 determines the new time 
averaged velocity of a characteristic droplet. However, whereas there is only one 
velocity per droplet, there is an entire temperature distribution per droplet. In the current 
analysis there are ten temperatures per droplet, corresponding to the ten nodes within the 
droplet. Therefore a time averaged temperature of each thermal node must be calculated. 
If T n in d is the temperature of the n th node of a thermal droplet entering a “pancake” and 
T n old d is the temperature of the n th node of a thermal droplet already in that “pancake”, 
then the new time averaged temperature of the characteristic droplet is: 


Eq. 56 


T n in d Ndrop in + (T n old d )Ndrop 
Ndrop jn +Ndrop 


[K] 


The SINDA VARIABLES 1 logic block of the droplet thennal model does a check to 
determine if there are any droplets residing in a given “pancake” at the current time step. 
As mentioned previously in the section, Model Details: Submodel “B”, Momentum , 
there may be moments in time that there are no droplets in a “pancake”. If there are no 
droplets in a “pancake”, the thennal nodes representing the droplet are placed in a 
boundary state. If droplets are present, the thermal nodes of the droplet are set to the 
distribution determined in the LLUINT LLOGIC 0 block determined by Equation 56. 
Also a heat rate from the condensation that occurs is applied to the outer most node of the 
droplet: 

Eq. 57 Q d = mcond h^ [W] 


It is not advisable to actually create a node that represents the condensate layer that forms 
around the droplet. This is because the layer is so thin and the thermal conductor so large 
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relative to the thermal conductors in the rest of the droplet that time step precipitously 
drops. 

If evaporation occurs, the heat rate on the droplet, becomes: 


Eq. 58 Q d = [ As s h conv (T*, - T satd ) - m eV a P h^ ] [W] 

If T sat is less than T satd : 

Eq. 59 Q d = [As s h C0nv (T 00 -T satd )-As s h filmd (T satd -T sat )(l-XL)-meva P h; g ] [W] 

Results 


Model results were compared to four Delta III upper stage hot fire tests that were run in 
the B2 facility. Figure 13 illustrates a summary table of the four hot fire runs. In all the 
cases presented below the droplets leaving the spray bar were 1500 microns in size and 
had an initial velocity 37 ft/sec. 



HOT 

HOT 

HOT 

HOT 


FIRE 3 

FIRE 6 

FIRE 8 

FIRE 10 

CONDENSING SPRAY CONDITIONS 





INLET CONDENSING SPRAY TEMPERATURE (DEG F)* 

50.6 

51.5 

55.99 

64.2 

INLET CONDENSING SPRAY FLOW RATE (KG/SEC) 

13878 

13878 

13878 

13878 

WATER LEVEL (FT) 

67.8 

73.8 

73.6 

64.5 

ULLAGE LENGTH (FT) 

45.65 

45.65 

45.65 

45.65 

ROCKET CONDITIONS 

ROCKET EXIT AREA (IN 2 ) 

1500 

1500 

1500 

1500 

ROCKET AREA RATIO 

77 

77 

77 

77 

ROCKET O/F RATIO 

6 

6 

6 

6 

ROCKET COMBUSTION PRESSURE (PSI) 

640 

640 

640 

640 


* For spray bar temperature rise due to engine heat exhaust or ejector heat output this was only an 
initial condition. 


Figure 13: Summary Table of Delta III Upper Stage Hot Fire Tests 

Figure 14 through Figure 17 illustrate individual test data cases compared to the respective 
SINDA/FLUINT model results. The models were first run with the following 
assumptions (i.e., Model ‘Hot Fire Case #’): 

• Steam formed from quench cooling along the rocket exhaust duct outside walls 
was neglected 

• Spray bar temperature rise from the heat of the engine exhaust. 

New model results were generated that included these two effects, (refer to output: 

Model ‘Hot Fire Case #’ New). 
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SHAY CKAUBtR PRESSURE VS TM 



Figure 14: Spray Chamber Pressure: Hotfire Test 3 and SINDA/FLUINT Model Results 
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Figure 15: Spray Chamber Pressure: Hotfire Test 6 and SINDA/FLUINT Model Results 
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Figure 16: Spray Chamber Pressure: Hotfire Test 8 and SINDA/FLUINT Model Results 


SPNAV CM* PRt$UJM£ VS 1M( 



Figure 17: Spray Chamber Pressure: Hotfire Test 10 and SINDA/FLUINT Model Results 


To examine effects of a candidate test article larger than the previously conducted engine 
tests at B-2, the SINDA/FLUINT Model was run using a two point engine test sequence 
lasting for 700 seconds. Figure 18 summarizes the significant input conditions used for 
this simulation. This model assumed the droplets leaving the spray bar were 1 500 
microns in size and had an initial velocity 37 ft/sec. The assumption that the spray bar 
water temperature rose due to the effect of engine exhaust heat was included in this 
analysis. Also vapor quench from cooling of the rocket exhaust duct was included. 
Figure 19 and Figure 20 provide the model outputs for spray chamber pressure and 
temperature. The useful result shown in these graphs are that the exhaust system , as 
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modeled, can support a 700 second duration engine firing where the initial engine thrust 
is higher followed by a lower thrust setpoint. The model can be useful in predicting 
exhaust system perfonnance for various hydrogen-oxygen engine combinations and 
testing durations. Future engine testing at B-2 will provide opportunities to evaluate and 
refine the model. 



Candidate Test 
Article, 

First 400 sec. 

Candidate Test 
Article, 

Last 300 sec. 

CONDENSING SPRAY CONDITIONS 

INLET CONDENSING SPRAY TEMPERATURE (DEG F)* 

40 

40 

INLET CONDENSING SPRAY FLOW RATE (KG/SEC) 

13878 

13878 

WATER LEVEL (FT) 

70 

70 

ULLAGE LENGTH (FT) 

49.25 

49.25 

ROCKET CONDITIONS 

ROCKET EXIT AREA (IN 2 ) 

5627 

5627 

ROCKET AREA RATIO 

243 

243 

ROCKET O/F RATIO 

5.797 

5.826 

ROCKET COMBUSTION PRESSURE (PSI) 

882 

637 


* For spray bar temperature rise due to engine heat exhaust or ejector heat output this was only an 
initial condition. 


Figure 18: Summary Table of Candidate Test Article Input 


SPRAY CHAMBER PRESSURE VS TIME 


Thrust 40466 Ibf 


Thrust 29225 Ibf 


TIME (SEC) 

Figure 19: Spray Chamber Pressure for Candidate Test Article SINDA/FLUINT 

Model Results 
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TEMPERATURE RISE (DEG R) 


CHAMBER SPRAY WATER TEMPERATURE RISE VS TIME 


TIME (SEC) 

Figure 20: Chamber Spray Water Temperature Rise for Candidate Test Article 

SINDA/FLUINT Model Results 
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Appendix A: CEA and SINDA/FLUINT Enthalpy and Entropy 

Reference States 


The relationship between the H 2 O enthalpies of REFPROP and CEA is as follows: 

Eq. A.l HcEA - Hf + (HreFPROP - HreFPROP satvap 298) [KJ/kg] 

where Hf is the enthalpy of formation of H 2 O vapor at 298. 15 K and 0. 1 MPA. It should 
be noted that the default setting for the reference state was set in REFPROP. Example 
A.l illustrates the use of Eq. 1 with the aid of Tables A.l through A. 3. 

Example A.l: Steam at 629.0 K and 0.1 MPA 

Hcea =-13423.647 + (3187.8 -2546.5) = -12782.3 [KJ/kg] 

Inspecting Table A.2, CEA’s output enthalpy of -12781.3 KJ/kg compares well to 
Example A. 1 . 

The relationship between the H 2 O entropies of REFPROP and CEA is as follows: 

Eq. A.2 ScEA - Sf + (SreFPROP - SreFPROP _satliq_298) [KJ/kg/K] 

where Sf is the absolute entropy of H 2 O liquid at 298. 15 K and 0. 1 MPA. Example A.2 
illustrates the use of Eq. 2 with the aid of Tables A.l through A. 3. 


Example A.2: Steam at 629.0 K and 0.1 MPA 

Scea = 3.882 + (8.4057 -0.36722) = 11.923 [KJ/kg/K] 


Inspecting Table A.2, CEA’s output entropy of 1 1.922 KJ/kg/K compares well to 
Example A.2. 


Table A.l: REFPROP Output of H 2 O Properties [Ref 5] 


TEMPERATURE 
K (SATURATED) 

PRESSURE 

MPA 

LIQUID 

DENSITY 

kg/m 3 

VAPOR 

DENSITY 

kg/m 3 

LIQUID 
ENTHALPY 
KJ/kg K 

VAPOR 

ENTHALPY 

KJ/kg 

LIQUID 

ENTROPY 

KJ/kg/K 

VAPOR 

ENTROPY 

KJ/kg/K 

298.15 

0.0031699 

997.0 

0.023075 

104.83 

2546.5 

0.36722 

8.5566 

TEMPERATURE 

K 

PRESSURE 

MPA 

DENSITY 

kg/m 3 

ENTHALPY 

KJ/kg 

ENTROPY 

KJ/kg/K 




629.00 

0.1 

0.34506 

3187.8 

8.4057 





38 


Table A.2: CEA Output of H 2 O Properties [Ref 6] 

THERMODYNAMIC EQUILIBRIUM PROPERTIES AT ASSIGNED 
TEMPERATURE AND PRESSURE 


STATE TEMP 
DEG K 

0.000 G 629.00 

PHI = 0.0000 

THERMODYNAMIC PROPERTIES 

P, MPA 0.10002 

T, DEG K 629.00 

RHO, KG/CU M 3.4455-1 

H, KJ/KG -12781.3 

U, KJ/KG -13071.6 

G, KJ/KG -20280.2 

S,KJ/(KG)(K) 11.9220 

M, MOLWT 18.015 

(DLV/DLP)T -1.00000 

(DLV/DLT)P 1 .0000 

CP, KJ/(KG)(K) 2.0348 

GAMMA (S) 1.2934 

SON VEL, M/SEC 612.7 

TRANSPORT PROPERTIES (GASES ONLY) 

CONDUCTIVITY IN UNITS OF MILLIWATTS/(CM)(K) 

VISC,MILLIPOISE 0.22572 

WITH EQUILIBRIUM REACTIONS 

CP, KJ/(KG)(K) 2.0348 

CONDUCITVITY 0.4969 

PRANDTL NUMBER 0.9243 

WITH FROZEN REACTIONS 


WT FRACTION ENERGY 
CHEMICAL FORMULA KJ/KG-MOL 

FUEL H 2.00000 O 1.00000 1.000000 

O/F = 0.0000 PERCENT FUEL = 1 00.0000 EQUIVALENCE RATIO = 1 .0000 


CP, KJ/(KG)(K) 2.0348 

CONDUCITVITY 0.4969 

PRANDTL NUMBER 0.9243 

MOLE FRACTIONS 

H20 1.00000 
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Table A.3: CEA Output of H 2 O Liquid Properties [Ref 6] 


THERMODYNAMIC EQUILIBRIUM PROPERTIES AT ASSIGNED 
TEMPERATURE AND PRESSURE 


CHEMICAL FORMULA 

FUEL H 2.00000 O 1.00000 

O/F = 0.0000 PERCENT FUEL = 100.0000 


WT FRACTION ENERGY STATE TEMP 

KJ/KG-MOL DEG K 

1.000000 0.000 L 298.15 

EQUIVALENCE RATIO = 1 .0000 PHI = 0.0000 


THERMODYNAMIC PROPERTIES 


P, MPA 

T, DEG K 
RHO, KG/CU M 
H, KJ/KG 

U, KJ/KG 
G, KJ/KG 

S, KJ/(KG)(K) 


0.10002 
298.15 
1.0799 6 
-15866.2 
-15866.2 
-17024.1 
3.8835 


M, MOL WT 
(DLV/DLP)T 
(DLV/DLT)P 
CP, KJ/(KG)(K) 
GAMMA (S) 

SON VEL, M/SEC 


-0.19929 

1.1114 

4.1821 

5.0179 

0.7 


MOLE FRACTIONS 


H20(L) 


1.00000 


The relationship between the H2 enthalpies of REFPROP and CEA, using the default 
reference setting, is not as straight forward. The CEA enthalpy is defined by: 

Eq. A .3 Hcea - (Hrefprop - Hassigned) [KJ/kg] 

noting that the heat of formation is zero. The assigned enthalpy, Hassigned, for liquid or 
vapor H2 is 9.012 KJ/mole or 4470 KJ/kg. The assigned enthalpy for superheated H2 is 
the enthalpy of H2 at 298.15 K and 0.1 MPA. However, for the fonner case, Hassigned 
should be “approximately” equal to the enthalpy of H2 at 298.15 K and 0.1 MPA plus the 
liquid enthalpy of H2 at 0. 1 MPA. The default reference setting of REFPROP does not 
yield this condition because the liquid enthalpy of H2 at 0.1 MPA is “almost” zero, (see 
Table 4). A new relative enthalpy value, at 298.15 K and 0.1 MPA, needs to be input 
(see Table 5). A new reference state for entropy, at 298.15 K and 0.1 MPA, also needs to 
be input. This value coincides with the value output by CEA (see Table 6 ). 


40 


Table A.4: REFPROP Output of H 2 Properties [Ref 5] 
Default Reference Setting 


LIQUID VAPOR LIQUID VAPOR LIQUID VAPOR 


TEMPERATURE 
K (SATURATED) 

PRESSURE 

MPA 

DENSITY 

kg/m 3 

DENSITY 

kg/m 3 

ENTHALPY 
KJ/kg K 

ENTHALPY 

KJ/kg 

ENTROPY 

KJ/kg/K 

ENTROPY 

KJ/kg/K 

20.27 

0.1 

70.805 

1.3368 

-0.073542 

445.40 

-0.0034731 

21.974 

TEMPERATURE 

PRESSURE 

DENSITY 

ENTHALPY 

ENTROPY 




K 

MPA 

kg/m 3 

KJ/kg 

KJ/kg/K 




298.15 

0.1 

0.081274 

3929.6 

53.436 





Table A.5: REFPROP Output of FI 2 Properties [Ref 5] 
Enthalpy Reference State = 4201 KJ/kg at 298.15 K and 0.1 MPA 


TEMPERATURE 
K (SATURATED) 

PRESSURE 

MPA 

LIQUID 

DENSITY 

kg/m 3 

VAPOR 

DENSITY 

kg/m 3 

LIQUID 
ENTHALPY 
KJ/kg K 

VAPOR 

ENTHALPY 

KJ/kg 

LIQUID 

ENTROPY 

KJ/kg/K 

VAPOR 

ENTROPY 

KJ/kg/K 

20.27 

0.1 

70.8 

1.337 

271.2 

716.7 

11.39 

33.36 

TEMPERATURE 

K 

PRESSURE 

MPA 

DENSITY 

kg/m 3 

ENTHALPY 

KJ/kg 

ENTROPY 

KJ/kg/K 




298.15 

0.1 

0.08127 

4201.0 

64.83 
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Table A.6: CEA Output of H 2 Properties [Ref 5] 

THERMODYNAMIC EQUILIBRIUM PROPERTIES AT ASSIGNED 
TEMPERATURE AND PRESSURE 


STATE TEMP 
DEG K 

0.000 G 298.15 

PHI = 0.0000 

THERMODYNAMIC PROPERTIES 

P, MPA 0.10002 

T, DEG 298.15 

RHO, KG/CUM 8.1337-2 
H, KJ/KG 1.1794 

U, KJ/KG -1228.54 

G, KJ/KG -19327.1 

S, KJ/(KG)(K) 64.8272 

M, MOLWT 2.016 

(DLV/DLP)T -1.00000 

(DLV/DLT)P 1.0000 

CP, KJ/(KG)(K) 14.3091 

GAMMA (S) 1.4050 

SON VEL, M/SEC 1314.4 

TRANSPORT PROPERTIES (GASES ONLY) 

CONDUCTIVITY IN UNITS OF MILLIWATTS/(CM)(K) 

VISC,MILLIPOISE 0.08915 

WITH EQUILIBRIUM REACTIONS 

CP, KJ/(KG)(K) 14.3097 

CONDUCITVITY 1.8274 
PRANDTL NUMBER 0.6981 

WITH FROZEN REACTIONS 

CP, KJ/(KG)(K) 14.3097 

CONDUCITVITY 1.8274 
PRANDTL NUMBER 0.6981 

MOLE FRACTIONS 

H2 1.00000 


WT FRACTION ENERGY 
CHEMICAL FORMULA KJ/KG-MOL 

FUEL H 2.00000 1.000000 

O/F = 0.0000 PERCENT FUEL = 1 00.0000 EQUIVALENCE RATIO = 1 .0000 
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Appendix B: Calculation of Droplet Inlet Velocity from the 

Condenser Spray Bar 


For an RL-10 engine test run in the B2 facility, the total flow rate injected through the 
condenser spray bar is: 

Eq. B.l msB =220000 [gal/min] = 13878 [lcg/s] 

The number of spray nozzles is: 

Eq. B.2 Nnoz = 694 

Therefore the flow rate per nozzle is: 

Eq. B.2 mnoz = 19.997 [kg/s] 

The nozzle exit radius is: 

Eq. B.4 Rnoz = 0.077458 [ft] = 0.02361 [m] 

Therefore the cross sectional area of a nozzle is: 

Eq. B.5 Acnoz= n(0.02361) 2 = 0.00175 [m 2 ] 

The velocity of the droplet exiting the nozzle is: 

• 

Eq. B.6 Vnoz= mnoz = 1 1 .427 [m/s] = 37.5 [ft/s] 

PsB ^Aioz 
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Appendix C: Abbreviations and Acronyms 


ANSI 

ASCI 

B2 

CEA 

CFD 

FORTRAN 

ft 

ft/sec 

gal/min 

H 2 

H 2 0 

ID 

J/kg/K 

K 

kg/m 3 

kg/s 

KJ/kg 

1/m/kg 

m 

m 2 

m/s 

MPA 

MSC Patran 


Spacecraft Propulsion Research Facility 
Chemical Equilibrium with Applications code 
Computational Fluid Dynamics 

feet 

feet per second 
gallons per minute 
hydrogen 
water 

Inner Diameter 

Joules per kilogram per degree Kelvin 
Kelvin 

kilograms per cubic meter 
kilograms per second 
Kilojoules per kilogram 
liters per meter per kilogram 
meters 

square meters 
meters per second 
Mega Pascal 


NASA National Aeronautics and Space Administration 

NIST 

OD Outer Diameter 

Pa Pascal 

REFPROP Reference Fluid Thermodynamic and Transport Properties 

sec seconds 

SMD Sauter Mean Diameter 

SINDA/FLUINT FLUINT is the fluid systems analyzer, SINDA is the thennal analyzer 
W watts 
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Appendix D: Obtaining and Running the Model 


SINDA/FLUINT 5.3 is very important. There are some double precision conversion 
issues at the moment. 
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Appendix E: Nomenclature 


k 

Thermal conductivity 

[W/m/K] 

g 

Gravity 

[m/s 2 ] 

h 

Heat transfer coefficient 

[W/m 2 /K] 

hconv HTUDIF 

Heat transfer coefficient for SINDA/FLUINT 
subroutine HTUDIF 

[W/nr/K] 

hfg 

Heat of vaporization 

[J/kg/K] 

h e ff 

Effective heat transfer coefficient between 
droplet surface and rocket exhaust duct flow 
saturation condition 

[W/m 2 /K] 

h f g 

Modified heat of vaporization to account 
for subcooling 

[J/kg/K} 

hfilm 

Film boiling heat transfer coefficient 

[W/nr/K] 

hquench 

Total quench heat transfer coefficient 

[W/m 2 /K] 

j 

Diffusion mass flux rate 

[kg/s/m 2 ] 

m 

Mass 

[kg] 


• 

m 

Mass flow rate 

[kg/s] 

• 

m water 

Quench water mass flow rate 

[kg/s] 

m 

Mass flux rate 

[kg/s/m 2 ] 

r 

Droplet radius 

[m] 

t 

Time 

[s] 

X 

Mass fraction 


Ac 

Cross-sectional area 

[m 2 ] 

As 

Surface area 

[m 2 ] 

CD 

Drag coefficient 


CF 

Correction factor for evaporation mass flow rate 


Cp 

Specific heat 

[J/kg/K] 

D 

Diameter 

[m] 

Do 

Outer diameter 

[m] 

D w h 

Diffusion coefficient of water and hydrogen 

[m 2 /s] 

F 

Friction factor 


FC 

SINDA/FLUINT FC factor 

[1/m/kg] 

FD 

Drag force 


Fg 

Body force 


G 

Thermal conductor 

[W/K] 

H 

Enthalpy 

[W/K] 

Hduct 

Enthalpy of shocked rocket exhaust 
duct flow 

[J/kg] 

Hquench 

Enthalpy of shocked rocket exhaust 

[J/kg] 
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Hwater 

duct flow with quench cooling water 
Enthalpy of quench cooling water 

[J/kg] 

Hcea 

duct flow 

Enthalpy output from CEA 

[KJ/kg] 

Hrefprop 

Enthalpy output from REFPROP 

[KJ/kg] 

HreFPROP satvap 298 

Saturated vapor enthalpy at 298.15 K 

[KJ/kg] 

Ja 

output from REFPROP 
Jakob number 


L 

Length 

[m] 

Nnoz 

Number of nozzles in spray bar 


Ndrop 

Number of droplets in a stratified tank 


Ndropi„ 

Number of droplets entering a stratified tank 


Nu 

Nusselt Number 


Ptot 

Total Pressure 

[Pa] 

Pr 

Prandtl number 


Q 

Heat rate 

[W] 

QL 

Heat rate applied to SINDA/FLUINT tank 

[W] 

R 

Radius 

[m] 

Re 

Reynolds number 


ScEA 

Entropy output from CEA 

[KJ/kg/K] 

Srefprop 

Entropy output from REFPROP 

[KJ/kg/K] 

SreFPROP satliq_298 

Saturated liquid entropy at 298. 15 K 

[KJ/kg/K] 

To 

Outer wall temperature 

[K] 

Tstat 

Static temperature 

[K] 

Tree 

Recovery temperature 

[K] 

Tref 

Reference temperature 

[K] 

Tstag 

Stagnation temperature 

[K] 

Twall 

Wall temperature 

[K] 

'-pn 

1 d 

Temperature of the n th node of a characteristic 

[K] 

T n in d 

droplet 

Temperature of the n th node of a characteristic 

[K] 

T n old d 

droplet entering a stratified tank 
Temperature of the n th node of a characteristic 

[K] 

V 

droplet in a stratified tank at the previous time 
step 

Velocity 

[m/s] 

Vind 

Velocity of characteristic a droplet 

[m/s] 

Voldd 

entering a stratified tank 
Velocity of characteristic a droplet 

[m/s] 

Vrel 

in a stratified tank at the previous time step 
Relative Velocity 

[m/s] 

XL 

Quality of a stratified tank 


P 

Density 

[kg/m 3 ] 
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Viscosity 


[kg/m/s] 


P 

Subscripts 

cond 

conv 

duct 

d 

evap 

h 

i 


n 

noz 

s 

sat 

satd 

satv 

satoo 

v 

walld 

w 

D 

SB 

i 

oo 


Condensation 

Sensible heat transport and transfer via 

convection and conduction 

Rocket exhaust duct 

Droplet 

Evaporation 

Hydrogen 

Interface between droplet condensate layer 
and diffusion layer (saturation conditions of 
droplet) 

Droplet thermal node number 
Spray bar nozzle 
Droplet surface 
Saturation 

Droplet saturation condition 
Saturated vapor 

Rocket exhaust duct flow saturation condition 

Vapor 
Droplet wall 
Water 

Duct 
Spray bar 

Liquid 

Exhaust duct flow 


Constants 


g 

n 

9.81 

3.1415 

Gravity 

[m/s 2 ] 

Hassigned 

4470.0 

Assigned enthalpy of liquid/vapor H 2 
in CEA 

[KJ/kg] 

H f 

-13423.647 

Enthalpy of formation of H?0 vapor 
at 298.15 K and 0.1 MPA 

[KJ/kg] 

s f 

3.882 

Absolute Entropy of H 2 O liquid 
at 298.15 K and 0.1 MPA 

[KJ/kg/K] 

P SB 

1000.0 

Liquid water density 

[kg/m 3 ] 
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